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1.  Introduction 


We  are  interested  in  simulating  a  variety  of  problems  in  3  dimensions  (3D)  fea¬ 
turing  large  electric  currents.  While  2D  simulations  have  been  quite  informative, 
cylindrical  symmetry  may  interfere  with  a  problem’s  relevant  physics.  Specifically, 
all  objects  in  the  domain  behave  as  if  they  are  extruded  360° — turning  particles  into 
hoops.  In  dealing  with  electrical  current,  this  can  have  serious  ramifications  on  the 
current  pathways.  In  3D  (r,  0,  z)  currents  can  adjust  their  pathways  anywhere  along 
those  360  degrees  given  the  right  conditions;  however,  in  2D  (r,  z)  those  pathways 
can  be  completely  choked  off  because  an  insulating  hoop,  rather  than  a  particle,  is 
present. 

Even  though  our  radial  trisection  mesh  is  a  90°  wedge  with  just  a  few  million  cells, 
in  developing  3D  exploding  wire  benchmarks,  we  found  that  simulations  were  pro¬ 
hibitively  expensive.  In  some  cases,  it  would  take  a  day  to  advance  1  ps.  Given  that 
we  require  up  to  approximately  100  ps  per  calculation,  this  is  clearly  a  problem. 

Users  have  the  ability  to  modify  a  large  number  of  solver  parameters,  some  of  which 
can  drastically  impact  performance.  However,  there  is  only  moderate  guidance  on 
what  parameters  to  tune  and  what  they  should  be.  Further,  it  is  expected  that  they 
can  be  problem  dependent.  Our  goal  then  is  to  quantify  these  uncertainties  and 
manually  optimize  the  problem  such  that  we  minimize  total  CPU  time  without  im¬ 
pacting  the  physics.  Only  recently  did  this  type  of  algorithmic  uncertainty  quantifi¬ 
cation  work  become  possible  through  a  competed  grant  of  high-priority  (Frontier) 
computing  time  from  the  DOD  High  Performance  Computing  Modernization  Pro¬ 
gram  supporting  this  and  other  research  areas. 

The  report  is  organized  as  follows.  Section  2  describes  the  computational  approach 
and  problem  setup  while  Section  3  provides  the  results  and  analysis.  Section  4  sum¬ 
marizes  the  key  findings  and  provides  recommendations. 

2.  Computational  Approach 

We  use  several  codes  and  algorithms  developed  by  Sandia  National  Faboratories 
(SNF)  to  evaluate  this  research.  AFEGRA  is  used  to  model  exploding  wires  via 
finite  elements  using  a  90°  wedge  cylindrical  mesh.  Simulations  were  run  on  Ex- 
calibur,  a  Cray  XC40,  and  for  comparison,  Haise,  an  IBM  iDataPlex.  We  found 
no  significant  difference  in  the  results.  AFEGRA  is  a  finite-element  multiphysics 
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simulation  code  developed  by  SNL  for  modeling  high-deformation  solid  dynam¬ 
ics,  shock-hydrodynamics,  magnetohydrodynamics  (MHD),  and  fcrroclcctricity. 1  1 
It  has  been  used  extensively  for  simulating  the  large  electric  currents  and  the  asso¬ 
ciated  heating,  material  effects,  and  MHD  phenomenology  of  exploding  wires.4-7 
The  MHD  module  of  ALEGRA  is  used  here,  which  incorporates  Maxwell’s  equa¬ 
tions  into  the  equations  of  solid  dynamics  using  the  “magneto-quasistatic”  approx¬ 
imation.8  In  this  approximation,  electromagnetic  radiation  and  displacement  cur¬ 
rents  are  neglected.  Further,  although  the  medium  is  assumed  to  have  resistivity  via 
Ohm’s  law,  and  does  carry  both  electric  current  and  electric  and  magnetic  fields, 
it  is  also  assumed  to  be  electrically  neutral — that  is,  the  electric  charge  density  is 
zero  everywhere,  and  the  medium  can  have  no  capacitive  behavior.  This  approxi¬ 
mation  is  sufficient  for  exploding  wires  in  the  regimes  considered  here,  so  long  as 
the  underlying  material  models  incorporate  the  effects  of  ionization  on  the  material 
electrical  conductivity.  For  exploding  wires,  the  Fee-More-Desjarlais  models  pro¬ 
vide  this  capability  and  is  known  to  be  reliable  through  the  solid,  liquid,  gaseous, 
and  warm  plasma  phases  for  common  metals.910 

In  cases  where  magnetic  forces  can  be  neglected,  one  may  use  the  “low  RM”  option, 
which  can  significantly  reduce  simulation  time.311  The  Magnetic  Reynold’s  num¬ 
ber,  Rm,  is  derived  from  the  equations  of  MHD  and  used  to  measure  the  relative 
importance  of  competing  physical  mechanisms.  It  is  defined  as  RM  =  HqctL2 /t, 
where  /io  is  the  permeability  of  free  space,  a  is  the  electrical  conductivity,  and  L 
and  r  are  the  characteristic  length  and  time  scales,  respectively,  of  the  system.  In 
our  case,  r(f)  oc  10~6  —  10“5,  a(t)  oc  102  —  106,  L(t)  oc  1CT3  —  10_1,  where  sev¬ 
eral  values  are  time-dependent.  At  early  times,  we  would  expect  RM  ~  1.  At  later 
times  when  the  wire  explodes,  this  simplified  analysis  is  complicated  by  the  fact 
that  conductivity  decreases  quickly  for  the  metal  plasma — not  for  the  remaining 
solid  wire  pieces — and  a  can  increase  for  the  medium.  Each  component  is  impor¬ 
tant.  Additionally,  the  length  scale  increases  as  the  plasma  expands,  and  this  term 
may  or  may  not  dominate  over  conductivity.  At  later  time,  therefore,  we  could  ex¬ 
pect  Rm  ~  0.1  —  10.  In  both  cases,  the  low  RM  option  (RM  «  1)  is  not  quite 
applicable  and  not  used  in  this  study. 

AFEGRA  uses  explicit  time  integration  to  solve  the  solid-dynamic  and  nonresistive 
component  of  the  multiphysics  equations,  which  are  hyperbolic.  However,  the  diffu¬ 
sion  of  magnetic  field  and  associated  resistive  heating  of  the  medium  are  captured 
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in  an  additional  system  that  is  elliptic.  These  equations  are  solved  implicitly  on 
each  timestep  using  iterative  methods  and  are  coupled  to  the  full  system  via  opera¬ 
tor  splitting.1  The  Trilinos/AztecOO  infrastructure  is  used  to  form  the  linear  system 
that  represents  this  problem  and  solve  it  using  a  variety  of  techniques  available  to 
the  user  via  input  file  controls. 

We  use  the  conjugate  gradient  (CG)  solver — an  iterative  method  for  solving  sparse 
linear  systems — as  implemented  in  Trilinos/AztecOO  with  a  user-specified  toler¬ 
ance  on  the  iterative  residual.  These  may  be  preconditioned  within  Trilinos/ML, 12,13 
and  since  we  are  solving  Maxwell’s  equations,  we  use  the  algebraic  multigrid  pre¬ 
conditioner.14  This  preconditioner  has  a  number  of  different  user-determined  set¬ 
tings  and  parameters.  These  settings  can  have  a  large  impact  on  simulation  time, 
and  it  is  not  always  clear  what  values  they  should  have.  Accordingly,  optimizing 
6  of  these  settings  is  the  primary  focus  of  this  report,  for  which  we  use  Dakota, 
a  toolkit  that  provides  algorithms  for  optimization,  uncertainty  quantification,  and 
sensitivity  analysis.15  It  can  be  called  directly  from  within  ALEGRA,  making  it  a 
very  attractive  tool.  In  this  report  we  use  it  primarily  to  manage  the  parameter  study. 

The  range  of  length  scale  in  most  exploding  wire  problems  is  enormous.  A  typical 
wire  diameter  can  be  as  small  as  tens  or  hundreds  of  microns,  but  the  expanding 
plasma  and  shock  wave  resulting  from  the  wire  burst  can  extend  to  tens  of  cen¬ 
timeters.  The  elliptic  nature  of  the  magnetic  diffusion  problem  also  demands  large 
computational  domains  so  that  boundary  conditions  are  as  near  to  infinitely  distant 
as  possible.  This  results  in  simulations  with  many  millions  of  elements  or  more, 
and  thousands  or  tens  of  thousands  of  timesteps.  Parallel  computing  platforms  are 
needed,  and  to  ensure  the  large-scale  parallel  performance  of  these  linear  solvers, 
the  multilevel  (ML)  package  from  the  Trilinos  library  is  used  in  ALEGRA,  coupled 
to  the  AztecOO  solvers. 

In  this  research,  we  seek  to  establish  optimal  settings  for  the  ML  configuration  in 
ALEGRA  exploding  wire  simulations.  The  ML  package  uses  algebraic  multigrid 
preconditioning  to  ensure  parallel  scaling.  As  such,  it  is  a  critical  part  of  ALEGRA- 
MHD  input  for  determining  the  overall  performance  and  rate  of  throughput  of 
ALEGRA-MHD  simulations. 

The  algebraic  multigrid  algorithms  in  ML  approximate  the  solution  of  a  linear  sys¬ 
tem  on  a  fine  mesh  with  a  series  of  solutions  on  progressively  coarser  meshes. 
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Figure  1  shows  a  sample  computational  mesh,  with  2  coarse  meshes  generated  by 
multigrid.  Each  coarser  mesh  represents  some  set  of  progressively  “lower  energy” 
or  “smoother”  components  of  the  solution  on  the  fine  level  solution.  By  decompos¬ 
ing  the  solution  in  this  fashion,  we  can  use  a  basic  iterative  method,  such  as  Jacobi, 
Chebyshev,  or  Gauss-Seidel  to  rapidly  reduce  error  that  is  “high  energy”  on  the  in¬ 
dividual  meshes.  For  the  purpose  of  this  study,  we  will  use  the  Chebyshev  smoother 
on  all  levels. 


Fig.  1  Example  of  a  multigrid  hierarchy,  with  meshes  named  A,  B,  and  Z  from  left  to  right 


Our  study  will  focus  on  6  parameters,  3  of  which  are  related  to  the  smoothing  pro¬ 
cess  described  above,  one  of  which  is  related  to  the  number  of  levels,  and  the  final  2 
require  additional  explanation.  The  first  parameter  we  consider  is  FINE  SWEEPS, 
which  represents  the  number  of  smoother  sweeps  used  on  the  user-supplied  fine 
mesh  (mesh  A  in  Fig.  1).  The  second  parameter  is  INTERMEDIATE  SWEEPS, 
which  represents  the  number  of  smoother  sweeps  used  on  all  meshes  except  the 
fine  and  coarse  meshes  (all  meshes  except  A  and  Z  in  Fig.  1).  The  third  parameter 
is  COARSE  SWEEPS,  which  represents  the  number  of  smoother  sweeps  used  on  the 
coarsest  mesh  (mesh  Z  in  Fig.  1).  The  fourth  parameter,  MULTILEVEL  GRIDS,  is 
the  maximum  number  of  levels  to  use.  Regardless  of  the  maximum  number  of  lev¬ 
els,  MF  will  stop  the  coarsening  when  the  coarse  grid  reaches  a  certain  (sufficiently 
small)  size. 

The  fifth  parameter  relates  to  MF’s  ability  to  coarsen  selectively.  In  the  case  of 
stretched  meshes,  or  nonconstant  material  parameters,  the  multilevel  preconditioner 
can  be  more  effective  by  coarsening  in  particular  directions  and  not  in  others.  In 
MF,  we  do  this  via  the  AGGREGATE  THRESHOLD  parameter.  We  use  MF’s  auxil¬ 
iary  aggregation  capability,  which  means  that  we  “drop”  connections  in  the  matrix 
corresponding  to  nodes  that  are  sufficiently  far  away  from  each  other.13 

The  final  parameter  involves  data  migration  between  processors  via  the  REPART  I T I  ON 
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parameter.  This  parameter  specifies  the  minimum  number  of  matrix  rows  to  require 
for  each  core.  If  the  number  of  rows  per  core  drops  below  this  level,  then  cores  are 
idled  to  keep  the  overall  work-to-communication  ratio  for  the  non-idle  cores  high. 


A  full  factorial  of  the  following  ML  parameters  and  their  ranges  are  investigated: 
FINE  SWEEPS  {1,5,10},  INTERMEDIATE  SWEEPS  {1,5,10},  COARSE 
SWEEPS  {1,5,10},  REPARTITIONING  {100,500,1000,2000}, 
MULTILEVEL  GRIDS  {5,10,15,20,25},  AGGREGATE  THRESHOLD  (AGG) 
{0.001,  0.005,  0.01,  0.05}.  This  consists  of  2160  simulations  where 
each  case  employs  roughly  180,000  cells  and  33  cores  per  simulation,  yielding  ap¬ 
proximately  5400  cells  per  compute  core.  The  parametric  study  was  managed  by 
the  Dakota  implementation  in  ALEGRA,  and  all  cases  ran  to  completion. 

Each  simulation  is  represented  by  a  90°  radial  trisection  mesh  populated  by  2  alu¬ 
minum  plates  and  an  aluminum  wire  with  a  sinusoidal  surface  that  oscillates  over 
25%  of  its  radius.  The  wire  material  density  is  initialized  as  a  random  field,  sam¬ 
pled  from  a  uniform  probability  distribution  spanning  ±1%  about  a  median  value  of 
2704  kg/m3,  with  an  aggregate,  or  “grain”,  size  that  is  twice  the  radial  cell  size.  This 
initial  density  field  is  stochastic  but  does  not  vary  from  one  simulation  to  the  next. 
Materials  are  situated  in  a  dry  air  environment,  the  latter  of  which  has  a  temperature 
ceiling  and  floor  of  50,000  K  and  50  K,  respectively.  For  these  simulations,  we  did 
not  perform  rigorous  convergence  work  (i.e.,  mesh  resolution,  domain  size)  because 
the  goal  was  to  only  find  an  optimum  setting  for  the  parameters  investigated.  We 
believe  those  parameters  would  not  deviate  significantly  under  such  convergence 
analysis.  A  wire  exploding  (rather  than  just  heating  up)  is  necessary  in  the  simu¬ 
lation,  as  we  expected  the  mixed  solid-liquid-plasma  environment  to  greatly  affect 
performance.  We  chose  a  larger-diameter  wire  because  smaller  diameters  would  re¬ 
quire  many  smaller  cells  to  resolve  it.  Those  smaller  cells  would  further  reduce  the 
time  step,  rendering  longer  simulation  times  and  becoming  very  expensive  given 
the  scope  of  this  study. 

Figure  2  illustrates  a  sample  simulation  where  the  wire  is  represented  by  an  isovol¬ 
ume  of  wire  temperatures  in  the  solid  phase.  The  wire  is  colored  by  current  density 
from  a  peak  value  down  6  orders  of  magnitude.  Note  in  the  wire’s  cross  section  that 
its  core  has  a  smaller  value  because  current  has  not  had  time  to  diffuse  inside  fur¬ 
ther.  The  wire  explosion  process  has  just  begun  at  the  smallest  diameters  of  the  wire, 
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where  the  temperatures  have  exceeded  that  of  melt  and  those  parts  of  the  isovolume 
are  no  longer  visible.  This  is  also  apparent  on  the  clip  plane  forming  a  backdrop  that 
plots  density  in  grayscale  over  roughly  5  orders  of  magnitude.  In  time,  the  event  is 
demarcated  by  a  vertical  line  in  the  plots  of  current  and  resistance — note  the  minor 
inflection  in  the  current  trace  assisted  by  a  large  increase  in  the  resistance  since  the 
liquid  and  gas  phase  for  aluminum  has  a  lower  conductivity  than  as  a  solid. 


Fig.  2  Sample  simulation  with  onset  of  wire  explosion  at  5.4  ps 

In  these  MHD  simulations,  total  CPU  time  can  be  greatly  reduced  by  setting  a  semi- 
relativisitic  limiter,  maxfast,  to  50  km/s.  This  limits  the  maximum  velocity  of  the 
fast  MHD  wave  speed,  which  is  directly  correlated  to  the  time  step  through  the 
Courant  condition.  Physically,  this  is  a  reasonable  choice  since  the  plasma  expan¬ 
sion  can  have  radial  velocities  up  to  5  km/s.  Our  effective  value  for  the  speed  of  the 
light  is  therefore  much  greater  than  relevant  velocities. 

Other  important  syntax  that  improves  simulation  stability  and  runtime  is  the  inclu¬ 
sion  of  (maximum)  temperature  clips — which  currently  require  a  SESAME  equa¬ 
tion  of  state  (EOS).  Some  EOS  tables  allow  for  exceedingly  high  temperatures  (108 
K).  Usually  when  temperatures  in  the  domain  exceed  105  K,  the  time  step  drops 
precipitously.  Rather  than  allowing  this  or  an  errant  cell  in  an  extreme  state  from 
controlling  and  possibly  crashing  the  simulation,  we  can  limit  the  maximum  allow¬ 
able  temperature.  One  must  use  this  with  care,  however;  if  the  temperature  clip  is 
too  low,  relevant  physical  processes  may  be  missed.  For  example,  we  found  that  a 
clip  of  25,000  K  for  all  material  was  over-constrained  since  increasing  it  affected 
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the  current  trace.  Also,  when  viewing  temperatures  spatially,  most  of  the  conduct¬ 
ing  region  was  clipping  at  the  maximum  value.  For  this  study,  we  settled  on  a  clip  of 
75,000  K.  For  nonelectrified  simulations  where  temperature-induced  phase  changes 
are  not  expected,  a  high-clip  of  3000-5000  K  is  reasonable.  Temperature  clips  are 
chosen  rather  than  the  more  heavy-handed  discard  option  of  deleting  material  from 
the  domain — based  on  certain  triggers — and  replacing  it  with  void.  This  alternate 
approach  violates  mass  conservation.  It  also  creates  a  further  problem  with  electro¬ 
magnetics  since  void  is  an  insulator  and  its  sudden  emergence  in  large  quantity  can 
choke  off  important  current  pathways.  This  problem  is  more  likely  to  emerge  in  2D 
than  3D,  but  we  still  believe  that  temperature  clips  are  more  appropriate  here. 

It  was  expected  that  changing  solver  parameters  should  have  no  effect  on  the  solu¬ 
tion.  We  verified  this  (not  pictured)  for  several  hundred  random  cases  and  found  the 
results  (e.g.,  current,  inductance,  and  resistance)  to  be  essentially  identical. 

3.  Results  and  Discussion 

Figure  3  illustrates  the  results  from  all  2160  simulations.  Each  subplot  shows  the  to¬ 
tal  CPU  time  as  a  function  of  the  values  for  the  corresponding  ML  parameters.  Note 
that  spaces  between  horizontal  gridlines  are  equivalent  to  roughly  5.5  h.  Black  dots 
represent  the  average  for  that  ensemble  of  data.  The  fastest  simulation  completed 
in  66,178  s  (18.4  h)  and  the  slowest  in  156,386  s  (43.4  h).  Clearly  Fine  Sweeps 
has  the  largest  effect  on  simulation  time,  with  the  smallest  value  being  the  best 
performer.  Over  this  full  range  of  inputs,  the  remaining  parameters  appear  to  show 
little  performance  trend. 
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Fig.  3  Timing  results  from  all  simulations 


Additional  performance  trends  emerge,  however,  if  we  look  at  the  fastest  simula¬ 
tions.  Accordingly,  we  selected  a  subset  of  data — 405  cases — where  CPU  <  80  ks 
(Fig.  4).  Fine  Sweeps  is  absent  in  this  figure  since  it  is  unity  for  all  of  these  simu¬ 
lations.  Note  that  spaces  between  horizontal  gridlines  are  now  roughly  equivalent 
to  1.5  h.  Finer  structure  is  now  apparent  in  both  intermediate  sweeps  and  aggre¬ 
gate  threshold,  where  in  both  cases  the  lowest  value  improves  performance  when 
FINE  SWEEPS=1.  There  do  not  appear  to  be  any  further  significant  correlations 
to  performance  among  the  remaining  variables  in  this  study. 
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Fig.  4  Timing  results  from  simulations  clocked  at  less  than  80  ks 

For  a  selection  of  data  where  AGG=0 . 001,  we  can  look  at  how  CPU  is  affected 
by  the  average  number  of  solver  iterations  per  cycle  (ITERSP)  and  the  average 
solver  solution  time  per  cycle  (SOLTIME)  in  a  given  simulation  (Fig.  5).  These 
are  user-defined  variables  rather  than  normal  variable  output.  In  both  cases  they 
represent  the  average  value  over  the  required  iterations  to  converge  for  each  cycle. 
In  the  left  panel,  we  see  a  clear  correlation  between  the  average  solution  time  and 
total  simulation  time.  However,  in  the  right  panel,  islands  of  data  document  a  sur¬ 
prising  result.  As  the  number  of  iterations  increases,  the  total  simulation  tends  to 
get  smaller.  This  is  counterintuitive  as  more  iterations  might  be  expected  to  require 
more  overall  time.  It  implies  that  the  cost  per  iteration  can  be  quite  low. 
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Fig.  5  Effects  of  average  solver  solution  time  and  iterations  on  total  solution  time 


When  using  the  Aztec  linear  solver,  one  must  always  find  a  balance  between  long 
and  costly  but  clever  iterations  (more  sweeps,  more  levels,  etc.)  versus  quicker, 
cheaper,  but  less  sophisticated  iterations  (fewer  sweeps,  fewer  levels).  For  some 
problems  the  former  gives  faster  overall  turnaround,  and  for  others,  the  latter  gives 
faster  turnaround.  Either  way,  the  problem  would  be  intractable  without  ML;  still, 
a  balance  has  to  be  found  to  minimize  the  time  to  solution.  The  data  on  the  right- 
hand  side  of  Fig.  5  suggest  that  for  the  3D  exploding  wire,  the  latter  option  is  more 
beneficial:  quicker,  cheaper  iterations.  It  is  unclear  what  the  conditions  are  for  the 
balance  to  swing  to  the  other  side  (i.e.,  a  much  smaller  wire  diameter,  or  perhaps  a 
mesh  with  a  much  larger  array  of  element  dimensions  and  aspect  ratios). 

4.  Conclusion 

We  performed  a  6-parameter  study  to  find  the  optimum  settings  in  the  Trilinos/AztecOO 
solver  for  minimizing  total  CPU  time  in  3D  exploding  wire  simulations  using  a  90° 
radial  trisection  mesh  in  ALEGRA.  We  found  the  best  performance  to  be  most  prob¬ 
able  when  using  FINE  SWEEPS=1 ;  INTERMEDIATE  SWEEPS=1; 
AGGREGATE  THRESHOLDS  .001,  where  the  former  has  the  largest  effect  on 
simulation  time.  A  separate  simulation  running  with  aggregate  threshold  set  at  5.0e- 
4  showed  no  further  improvements.  As  expected,  changing  solver  parameters  had 
no  effect  on  the  solution.  We  verified  this  for  several  hundred  random  cases  and 
found  the  results  (e.g.,  current,  inductance,  and  resistance)  to  be  essentially  identi¬ 
cal.  Our  complete  input  then  for  the  solver  is  as  follows: 
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aztec  1 
solver,  eg 
scaling,  sym_diag 
conv  norm,  rhs 
tol=  l.e-6 
max  iter  =  1000 
output  =  last 
multilevel 

verbose,  1 

formulation  =  maxwell2 
fine  sweeps  =  1 
intermediate  sweeps  =  1 
coarse  sweeps  =  1 
multigrid  levels  =  25 

interpolation  algorithm  =  uc  aggregation 
repartition,  minimum  edges  per  processor  =  500 
smooth  prolongator 
aux  aggregation 
aggregate  threshold  0.001 

end 

end 

Tighter  tolerances  will  of  course  slow  down  performance  and  may  be  required  for 
certain  problems  to  properly  converge.  Turning  off  AUX  AGGREGATION  had  no 
effect  on  performance  time. 

Total  simulation  time  for  this  class  of  MHD  problems  can  be  further  improved 
with  the  following  additions.  Dominant  among  these  is  SEMIRELATIVISTIC, 
MAXFAST,  which  limits  the  fast  MHD  wave  speed.  We  found  that  setting  it  at  5- 
10  times  the  maximum  physically  relevant  velocities  (i.e.,  the  plasma  expansion 
velocity)  greatly  improved  performance.  However,  formal  convergence  analysis  is 
still  required  to  find  a  minimum  acceptable  value.  Additionally,  we  recommend 
that  high  clips  for  temperature  be  used  when  possible  for  added  stability  and  speed 
to  prevent  an  occasional  grid  cell  from  reaching  cripplingly  large  temperatures, 
usually  due  to  mixed  materials.  Discarding  material  is  also  an  option  once  extreme 
states  are  reached;  however,  it  introduces  void,  an  insulator,  which  can  alter  the 
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current  path  or  even  quench  it  if  enough  void  is  present.  In  these  simulations  we  set 
Tmax  =  75,  000  K.  A  density  floor  of  0.01  kg/m3  is  a  reasonable  value  and  helps 
prevent  the  time  step  from  dropping  too  low.  All  of  these  enhancements  resulted  in 
up  to  a  10-fold  improvement  in  total  CPU  time. 
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List  of  Symbols,  Abbreviations,  and  Acronyms 


1-2-3D  1 -2-3-dimensional 
CG  conjugate  gradient 

CPU  central  processing  unit 

EOS  equation  of  state 
MHD  magnetohydrodynamics 
ML  multilevel 

SNL  Sandia  National  Laboratories 
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